Dietary characterization of the endangered salt marsh harvest mouse and sympatric rodents using DNA metabarcoding

Abstract The salt marsh harvest mouse (Reithrodontomys raviventris; RERA) is an endangered species endemic to the coastal wetlands of the San Francisco Estuary, California. RERA are specialized to saline coastal wetlands, and their historical range has been severely impacted by landscape conversion and the introduction of non‐native plant and rodent species. A better understanding of their diet is needed to assess habitat quality, particularly in relation to potential competitors. We investigated three questions using DNA metabarcoding with ITS2 and trnL markers: (1) Do RERA specialize on the native plant, pickleweed (Salicornia pacifica), (2) Do RERA consume non‐native plants, and (3) What is the dietary niche breadth and overlap with three sympatric native and non‐native rodents? RERA diet was dominated by two plants, native Salicornia and non‐native salt bush (Atriplex spp.), but included 48 plant genera. RERA diet breadth was narrowest in fall, when they consumed the highest frequencies of Salicornia and Atriplex, and broadest in spring, when the frequencies of these two plants were lowest. Diet breadth was slightly lower for RERA than for co‐occurring species in pairwise comparisons. All four species consumed similarly high frequencies of wetland plants, but RERA consumed fewer grasses and upland plants, suggesting that it may be less suited to fragmented habitat than sympatric rodents. Diet overlap was lowest between RERA and the native California vole (Microtis californicus). In contrast, RERA diet overlapped substantially with the native western harvest mouse (R. megalotis) and non‐native house mouse (Mus musculus), suggesting potential for competition if these species become sufficiently abundant.


| INTRODUC TI ON
The salt marsh harvest mouse (Reithrodontomys raviventris; RERA; Figure 1) is a habitat specialist occurring solely in the salt marshes of the San Francisco Estuary (SFE), California, USA. They are the only known mammal species restricted to coastal marshes (Greenberg et al., 2006). Despite being listed as an Endangered Species since the inception of the US Endangered Species Act (United States Fish and Wildlife Service [USFWS], 1970), the ecology of RERA remains poorly known. Most prior research effort has emphasized habitat associations (Smith et al., 2018a). Other aspects of RERA ecology, including diet, predation, disease ecology, and interspecific interactions, remain poorly understood (Smith et al., 2018b). Historically, RERA has been considered a specialist of Salicornia marsh habitat (Fisler, 1965;Shellhammer et al., 1982;USFWS, 1984). Recent evidence, however, suggests that RERA may be less specialized to Salicornia habitat than believed, particularly in brackish marshes with lower salinity and greater plant diversity (Smith et al., 2020;Smith & Kelt, 2019;Sustaita et al., 2011). The SFE has been altered by over a century of anthropogenic impacts, including the loss of >90% of historical tidal marsh habitat (Hobbs et al., 2006;Williams & Faber, 2001). Reflecting these threats, RERA is listed as endangered by the state of California (California Department of Fish and Wildlife [CDFW], 1971), the federal government (USFWS, 1970), and the International Union for Conservation of Nature (IUCN, 2021).
Three other rodent species are commonly detected in SFE marshes, the western harvest mouse (R. megalotis; REME), California vole (Microtus californicus; MICA), and house mouse (Mus musculus; MUMU). REME is native to this region, ranges throughout the western United States, and is considered a habitat generalist (Webster & Jones, 1982). REME are less salt tolerant than RERA and are thought to occur primarily in uplands and marsh-upland edge in the SFE (Fisler, 1965). MICA also is native and is considered a grassland specialist (Batzli & Pitelka, 1971). MICA are notably larger than both harvest mouse species, and likely are behaviorally dominant to them, although RERA may be better adapted to high salinity conditions (Blaustein, 1980;Geissel et al., 1988). Previous work has indicated both negative (Geissel et al., 1988) and positive (Sustaita et al., 2011) associations between RERA and MICA habitat use, but their diet interactions remain unknown. Finally, nonnative MUMU commonly co-occur with RERA throughout the SFE (e.g., Bias & Morrison, 2006;Marcot et al., 2020). MUMU are highly fecund (Bronson, 1979;Pye, 1993), opportunistic, and tolerant of a wide range of ecological conditions (Berry, 1981). MUMU and RERA are compatible in captivity (Catlett & Shellhammer, 1962), and may (Bias & Morrison, 2006) or may not (Sustaita et al., 2011) partition habitat, but the recapture probabilities for MUMU were positively influenced by RERA densities at Suisun Marsh (Smith et al., 2020).
Non-native species can affect multiple trophic levels of a community, as they may represent novel predators, novel competitors, or novel food resources to different community members (Lepczyk & Rubinoff, 2017). Ecological specialists may be particularly sensitive to non-native species, as they may be less likely than generalists to utilize novel resources (Abernethy et al., 2016;Marvier et al., 2004).
Indeed, the ecological impacts of non-native species are considered the leading cause of extinction for endemic mammals (Pimm et al., 1995). Given the challenges that ecological specialists must overcome to persist in a highly altered ecosystem, understanding the effects of non-native food resources and non-native intraguild species are critical ecological underpinnings for RERA conservation.
Despite the importance of strong ecological baselines to conservation and management of endangered species, the dietary habits of RERA remain poorly understood. RERA diet has been inferred from habitat use (USFWS, 2010), characterized coarsely by stomach content analysis (Fisler, 1965), and measured with cafeteria trials (Smith & Kelt, 2019). Rather than providing consensus, these studies have led to divergent views of RERA diet. Habitat associations and stomach content analyses suggest that RERA consume primarily Salicornia (Fisler, 1965). Conversely, cafeteria trials suggest that RERA are generalist foragers that may prefer some non-native plant species over Salicornia. Additionally, comparative dietary interactions of RERA and sympatric rodents have never been investigated.
To address these critical knowledge gaps, we applied DNA metabarcoding to fecal samples collected from rodents in the SFE.
Metabarcoding often identifies significantly more dietary taxa at finer taxonomic levels than other methods (Kartzinel et al., 2015;Soininen et al., 2009;Valentini et al., 2009). Additionally, the noninvasive nature of dietary metabarcoding makes it particularly appealing for research on threatened and endangered species (e.g., Castle et al., 2020;Iwanowicz et al., 2016). Our objectives were to describe diets of RERA and sympatric rodents and characterize both spatial and temporal dietary variation. We evaluated the hypothesis that RERA has a more specialized diet than sympatric species, and assessed the potential for competition over food resources.
F I G U R E 1 Salt marsh harvest mouse (Reithrodontomys raviventris) in salt marsh habitat in Suisun marsh, California, USA 2 | ME THODS

| Sample collection
We collected feces from animals captured during regular livetrapping surveys (see Smith et al., 2020 for details of survey design and associated protocols). Preliminary trials showed that feces collected directly from live-trapped animals were more likely to be composed entirely of diet items from the trapping bait, so we endeavored to collect feces from the bedding in traps to characterize diet before consumption of bait. Samples were collected at five sites in coordination with the California Department of Fish and Wildlife during regular RERA monitoring ( Figure 2). One of these sites-the Goodyear Slough Unit (GYS) of the Grizzly Island Wildlife Area-was trapped quarterly over 2 years (Summer 2018-Spring 2020, inclusive), allowing us to partition diet into four seasonal data sets; all other sites were trapped once either in summer or late spring, resulting in a total of eight sampling units ( Table 1). We also surveyed vegetation plots to characterize availability of potential diet items at

| Laboratory procedures
We extracted DNA from fecal samples using Qiagen Plant Mini Kits (Qiagen, CA, USA). For each captured individual, we extracted DNA from pooled fecal pellets; we targeted >5 pellets from each individual, and final pellet numbers in extractions ranged from 1-13 (mean = 5.7). Library preparation followed the general template of the Illumina 16S metagenomic protocol (Illumina, 2015). Since single markers may only amplify a subset of plant taxa in herbivore diets (Goldberg et al., 2020), we applied two commonly used plant metabarcoding markers. We amplified the second internal transcribed spacer (ITS2), which is a longer fragment (~290-340 bp) of nuclear ribosomal DNA with high taxonomic resolution (China Plant Barcode of Life Group et al., 2011), and the P6 loop of the trnL intron, a shorter fragment (~25-90 bp) of chloroplast DNA, which is less likely to be affected by degradation but has coarser taxonomic resolution (Fahner et al., 2016). We used the R package "PrimerMiner" (Elbrecht & Leese, 2016) (Cheng et al., 2016) for ITS2, and the primers F I G U R E 2 Map of sampling locations in this study. Five sites were livetrapped for salt marsh harvest mice (Reithrodontomys raviventris) and other small mammals and fecal pellets were collected to characterize diet with DNA metabarcoding. Sites included Goodyear Slough (GYS) and crescent unit (CRES) of the Grizzly Island wildlife area, ponds 1&2 (HS12) and area 9 (HS9) of the Hill Slough wildlife area, and EDEN landing ecological reserve (EDEN). GYS was sampled quarterly over two years to provide seasonal dietary data. All other sites were trapped opportunistically on one occasion each, either in late spring (CRES and EDEN) or summer (HS12 and HS9) trnl_g and trnl_h (Taberlet et al., 2007) for trnL. We added sequence overhangs to the 5′ ends of amplicon primers to facilitate annealing to Illumina sequencing adapters (compete primer sequences in Appendix A1, PrimerMiner scores in Appendix A2). We amplified ITS2 using the thermal protocol described in Moorhouse-Gann et al. (2018) and amplified trnL using the thermal protocol described in Taberlet et al. (2007). Given that biological replication (i.e., samples from unique individuals) yields significantly more variation in diet than technical replication (i.e., multiple PCR replicates per individual), we chose to prioritize our resources for biological replication and therefore conducted a single PCR replicate for each individual (Mata et al., 2019).
We included 20 positive controls and 12 negative controls per sequencing lane (Appendix A3). Positive controls were composed of DNA extracted from plants collected from our field sites. Each set of 20 positive controls included 10 single-species controls to assess sensitivity and to help estimate misassignment error based on the proportion of nontarget reads within single-species controls. We also included 10 two-species controls, which had equal concentrations of DNA from two plant taxa and helped to understand potential amplification biases. We used deionized water for negative controls.
We sequenced libraries using MiSeq 300 PE for ITS2 and 75 PE for trnL. Sequencing and sample demultiplexing were conducted by the UC Davis Genome Center.

| Bioinformatic processing
We trimmed and quality-filtered sequences using cutadapt (Martin, 2011). We identified Amplicon Sequence Variants (ASVs) using DADA2 (Callahan et al., 2016). To identify the taxonomy of ASVs, we created a custom database of ITS2 and trnL sequences of all plant genera known to occur in the SFE using the batch_download feature of PrimerMiner, which obtains sequences from both NCBI (Benson et al., 2013) and BOLD (Ratnasingham & Hebert, 2007) databases, and we manually re-formatted the reference sequences for use in DADA2. We used the assignTaxonomy feature of DADA2 to assign ASVs against the custom database, and used BLAST (Zhang et al., 2000) to corroborate assignments. We assigned ASVs at the genus level, except for some trnL sequences that could not be assigned to a single genus and were therefore assigned to the lowest possible suprageneric level (e.g., family or multiple genera).
We conducted sequence processing and assignment independently for each MiSeq lane. After taxonomic assignment, we retained only ASVs that comprised >0.01% of the total sequence reads in a lane. We then used positive and negative controls to inform filtering parameters to account for misassigned ASVs (O'Rourke et al., 2020). Based on the negative and positive controls, we discarded any sample with <5000 (trnL) or <3000 (ITS2) sequencing reads, or with <20% of reads successfully matching plant taxa; and within samples, we discarded any taxa comprising <0.5% (trnL) or <1.0% (ITS2) of reads. After applying those filters, we removed any taxa that likely originated from a source other than wild RERA diet (Avena, Helianthus, Juglans, Panicum, and Phalaris from trapping bait, and Gossypium from trap bedding).

| Salt marsh harvest mouse diet
We recorded the presence/absence of diet taxa within individual diets. We calculated the frequency of occurrence (FO = the proportion of individuals that consumed a given diet item) of diet items within RERA samples pooled across all sampling locations and seasons. We then categorized diet items as native versus non-native, by life form (grass, shrub, forb, vine), and by habitat (e.g., whether they were typical of wetlands or of uplands; determined from Jepson eFlora [Jepson Flora Project, 2021]), and we estimated FO for each category. We chose to use the presence/absence-based data (i.e., FO) due to the complexity of estimating biomass from relative read abundance (e.g., Deagle et al., 2019) with two different markers that TA B L E 1 Number of individuals of four rodent species captured during each season at each site: Salt marsh harvest mouse (Reithrodontomys raviventris; RERA), western harvest mouse (R. megalotis; REME), house mouse (Mus musculus; MUMU), and Note: Goodyear Slough (GYS) was surveyed during all four seasons (Su = summer, fa = fall, Wi = winter, and Sp = spring), whereas Hill Slough 1&2 (HS12), Hill Slough 9 (HS9), crescent unit (CRES), and EDEN landing (EDEN) were sampled only in summer or late spring. Due to varying sample sizes by site and season, interspecific comparisons of diet were conducted on a pairwise basis between RERA and one other species independently, and only at sites where diet data were available for ≥4 individuals of both species. a Sites with sufficient sample size to be included in comparisons of RERA and REME diet. We evaluated seasonal variation in RERA diet at Goodyear Slough. For this analysis, we pooled data by season across the 2 years, and estimated FO of diet items within each season. We tested for significant seasonal differences using the anosim function in the "vegan" R package (Oksanen et al., 2020). We compared diet to plant availability with Manly's Selection Index (W i ;Manly et al., 2002) using the R package "adehabitatHS" (Calange, 2011). We considered diet items to be "selected" when W i ± 95% confidence intervals >1, and "avoided" when W i ± 95% CI <1. Since vegetation availability data were collected at the genus level, we excluded any diet items identified at a coarser taxonomic level from selection analyses. We quantified dietary niche breadth as the effective number of species ( 1 D; Hill, 1973;Chao et al., 2014) derived from Shannon's Diversity Index (Shannon & Weaver, 1949) for unequal sample sizes and presence/absence data (Chao et al., 2014). We used the R package "iNEXT" (Chao et al., 2014) to estimate 1 D and 95% CIs using 500 bootstrap replicates. We considered seasonal differences in dietary niche breadth significant if 95% CIs were non-overlapping.
Additionally, we estimated dietary niche overlap between pairs of seasons using Jaccard's Similarity Index (J s ), calculated in "vegan." To visualize dietary niches in ordination space, we conducted nonmetric multidimensional scaling (nMDS) based on Jaccard Distances (J D ) using the metaMDS function in "vegan," and calculated 95% confidence ellipses for each season. We conducted nMDS over a range of dimensions (k) and selected the minimum number of dimensions (k = 3) in which stress of the ordination was <0.10.
We evaluated spatial variation in RERA diet from five sites sampled in late spring and summer (henceforth, "summer" sampling

| Salt marsh harvest mouse diet
We documented 53 taxa, including 48 genera and 5 higher order and Salicornia (Table 2B).

| Comparison of diet to co-occurring rodents
The FOs of REME and MUMU were qualitatively similar to those of RERA ( Figure 5; RERA had significantly lower dietary niche breadth than all three sympatric rodents ( Figure 6). Finally, dietary niche overlap in nMDS space was very high with both REME and MUMU, which effectively subsumed RERA dietary niche space (Figure 7a,b). In contrast, ordination highlighted that the diets of RERA and MICA were effectively distinct ( Figure 7c).

Salt marsh harvest mouse diet was dominated by Salicornia and
Atriplex year-round, but also included a wide variety of other native and non-native plants. Seasonal niche breadth was narrowest in fall when they consumed primarily Salicornia, Atriplex, and Grindelia.
RERA diet was less diverse than the diets of sympatric rodents due to less frequent consumption of grasses and upland plants. Equally important, however, is that Atriplex was nearly as prominent in RERA diet as Salicornia. A primary difference between the two plants was the lower availability of Atriplex, which led to relatively less consumption overall but high selection coefficients. Overall, these data support a hypothesis that Salicornia stands including mixtures of plants such as Atriplex, Frankenia, and Grindelia may provide more value to RERA than those with Salicornia alone (Fisler, 1965;Shellhammer et al., 1982). In particular, a growing body of work from Suisun Marsh, where brackish water promotes more diverse plant communities, has emphasized the importance of mixed vegetative communities over Salicornia-dominated sites (Botti et al., 1986;Smith & Kelt, 2019;Sustaita et al., 2011). Our data clarified that Salicornia is an important element in RERA diet, but that their diets were not strictly specialized.

| Seasonal changes in salt marsh harvest mouse diet
Optimal foraging theory suggests that animals will specialize on preferred foods when they are available, and that they will broaden their diets when preferred foods are unavailable (MacArthur & Pianka, 1966;Stephens & Krebs, 1986). In fall, RERA diet narrowed sharply and was overwhelmingly composed of three species (Salicornia, Atriplex, and Grindelia). In spring, however, consumption of these three plants declined and their dietary breadth expanded accordingly. We suspect that RERA foraging patterns may largely be driven by affinities for these three plants.
Dietary seasonality, in turn, likely is driven by plant phenology.
Fall, when RERA diet narrowed to focus almost exclusively on Salicornia, Atriplex, and Grindelia, is the peak seeding period for these three plants (Hutchings & Russell, 1989; Jepson eFlora Project et al., 2021), and is followed by dormancy or dieback in late winter and early spring, which coincided with reduced consumption by RERA. Whereas annual dieback of Atriplex has led some to suggest that this plant has limited value to RERA in winter and spring (Botti et al., 1986;USFWS, 2010), our data suggest substantial consumption of Atriplex year-round despite seasonal dieback. In contrast, some non-native plants, such as Phragmites, were consumed primarily during one season, and were avoided TA B L E 2 Manly's selection index (W i ) for plant genera in salt marsh harvest mouse (Reithrodontomys raviventris) diets in (a) four seasons at Goodyear Slough, and (b) summer at 5 locations/sites: Goodyear Slough (GYS), Hill Slough 1&2 (HS12), Hill Slough 9 (HS9), crescent unit (CRES), and EDEN landing (EDEN)  (Johnston, 1957). We observed higher frequencies of Grindelia in RERA diets during fall and winter, which may reflect an increase in habitat use associated with seasonally high tides. Taken together, these observations suggest that Grindelia may provide an important combination of high tide refuge, cover from predators, and forage to RERA during extreme diurnal high tides of late fall and early winter.
The diet of RERA broadened in spring, with increased consumption of upland plants that were negligible in the diets in other seasons. This was particularly notable for upland grasses, which is consistent with previous RERA stomach content analyses (Fisler, 1965). RERA remain largely restricted to marsh habitat with the exception of spring forays into terrestrial grasslands (Shellhammer et al., 1982;USFWS, 2010;Zetterquist, 1977). Fisler (1965) speculated that vegetative cover in grasslands was insufficient for RERA outside the spring growing season. Geissel et al. (1988) suggested that RERA retreated to uplands in response to springtime population irruptions of larger bodied voles.
Although we cannot discern whether competition or seasonal resource exploitation drove this pattern, our data support the hypothesis that utilization of terrestrial grasslands by RERA is largely limited to spring.
Several seasonal patterns in our data mirrored observations from cafeteria trials (Smith & Kelt, 2019). Seasonal selection indices of Salicornia and Atriplex in our study were high in fall and low in spring, corresponding with seasonal preference rankings in cafeteria trials.
In contrast, our data showed high FO of these plants in summer as well, whereas feeding trials did not. Our data also aligned with feeding trials that suggested increased preference for annual grasses in spring. On the other hand, feeding trials suggested high or moderate preference for Juncus in multiple seasons; we found low consumption of Juncus both overall and in proportion to availability in the present study. Despite high availability and high FO in MICA diet, we detected Juncus in just one of 189 RERA samples at Goodyear Slough. Another major conclusion from feeding trials was a strong preference for Polypogon. Our ability to corroborate this result may have been limited by low availability or absence of this plant from most of our study sites. Polypogon was rare at Goodyear Slough and relatively common at Crescent Unit, and consumption by RERA occurred in proportion to its availability.

| Dietary comparisons to co-occurring rodents
Diet of the endangered RERA overlapped substantially with that of the widespread REME, driven primarily by high frequencies of Salicornia and Atriplex. Although the kidney physiology of REME suggests capability to consume Salicornia, they were unable to survive in feeding trials after consuming even small amounts of this plant   (Coulombe, 1970). Similarly, captive REME starved when presented with only Salicornia and Distichlis as food sources (Fisler, 1965).
Nonetheless, our data revealed that wild REME regularly consumed both of these genera ( Figure 5, Appendix A9). REME were the only species to consume Atriplex (which grows primarily in diked wetlands) more frequently than Salicornia (which occurs frequently in both diked and tidal wetlands). This pattern most likely reflects differential space use, as REME are more abundant on diked wetlands than tidal wetlands, and RERA and MUMU abundances do not differ among wetland types (Smith et al., 2020). REME also consumed grasses (Distichlis, Festuca, Hordeum, and Phragmites) and upland plants such as thistles (Sonchus and Cynareae) with greater frequency than did RERA.
We also documented considerable dietary overlap between MUMU and RERA, driven by high frequencies of Salicornia and Atriplex. Similar to REME, MUMU consumed more grasses (Distichlis, Phragmites, Hordeum, and Festuca) and upland plants (Sonchus, Cynareae) than did RERA. In studies of habitat use in the SFE, MUMU were more closely associated with terrestrial grasses and fragmented habitat assemblages than were RERA (Bias, 1994;Bias & Morrison, 2006). Interestingly, despite a relatively generalist diet, only a single house mouse (of 26) consumed Grindelia, which was one of the most frequently consumed plants for the three native rodents ( Figure 5, Appendix A9).
Relative to RERA, the most distinct diet was that of MICA, primarily due to reduced use of Atriplex and a high frequency of Juncus.  (Batzli & Pitelka, 1971). In fact, the diet of MICA in our study more closely resembled that of Amargosa voles (M. c. scirpensis), a subspecies endemic to wetlands in the Mohave Desert (Castle et al., 2020), than MICA occupying the uplands adjacent to SFE marshes (Batzli & Pitelka, 1971).
The diet of RERA was more restricted ( Figure 5) and significantly less diverse ( Figure 6) than that of sympatric species.
Preference for Salicornia and Atriplex was notably greater for RERA, while sympatric species consumed higher proportions of several other species ( Figure 5). Notably, many of these latter plants were restricted to uplands, indicating that sympatric rodents are better equipped than RERA to utilize resources in edge habitats. Indeed, we note that REME, MUMU, and MICA generally are considered upland species, thus, our characterizations of their diets is specific to the individuals occurring on the upland/marshland edges and likely not reflective of these species as a whole.
Habitat fragmentation and small patch size reduce the probability of RERA occurrence (Bias & Morrison, 2006;Marcot et al., 2020), and occupancy of marsh habitat by REME and MUMU may be dependent upon the degree of habitat fragmentation and penetration of terrestrial grass microhabitats into the marsh (Bias & Morrison, 2006;Fisler, 1965). Our results support these important management issues, adding to a growing literature suggesting that fragmentation of marsh habitat and the associated increase in edge habitat are potential threats to RERA with respect to competition from upland-adapted rodents.

| CON CLUS IONS
We characterized the diet of RERA and three sympatric rodents in remnant coastal marsh habitat of the SFE. Salicornia and Atriplex were prominent in RERA diet across sites and seasons. RERA diet narrowed sharply in fall during peak seed production of Salicornia, Atriplex, and Grindelia, which appeared to be favored food items.
RERA consumption of terrestrial grass was largely restricted to spring, coinciding with previously documented patterns of seasonal use of upland habitats. RERA diet overlapped substantially with REME and the non-native MUMU, but not with the native MICA.
Our data provide the first comprehensive characterization of RERA diet in the wild. This information fills critical knowledge gaps in the ecology of RERA and can guide habitat and vegetation management decisions to benefit conservation of the species. Moreover, our study lays the groundwork for future investigation of competition affecting this endangered species.

F I G U R E 7
Population-level dietary overlap as represented in ordination (nMDS) plots of (a) salt marsh harvest mouse (Reithrodontomys raviventris; RERA) compared with western harvest mouse (R. megalotis; REME), (b) RERA compared with house mouse (Mus musculus; MUMU), and (c) RERA compared with California voles (Microtus californicus; MICA). Dots represent population-level diet using frequency of occurrence data. Ellipses show 95% confidence intervals